
##create dataset that drops after and before treatments
limited<-subset(study3, treatment<4)
table(study3$treatment
      )
#Compare means of the four items by treatment
#Plot the Means
library(ggplot2)
library(Rmisc)
library(DescTools)
library(cowplot)
library(effectsize)

##make a better treatment indicator

limited$treatment2<-as.factor(limited$treatment)
levels_temp=c("Control", "Boy", "Girl", "Before", "After")




##Plot the means by condition
RES_mean <- summarySE(limited, measurevar="Respect", groupvars="treatment2",
                      na.rm=TRUE, conf.interval=.95)
RES_mean$test<-c("Control", "Boy", "Girl")
RES_mean


plot_respect<-ggplot(RES_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), y=Respect, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Respect-ci, ymax=Respect+ci), colour="red", data=RES_mean) +
  geom_point(shape=21, size=1, fill="white")+
  xlab("Treatment")+
  ylab("Proportion choosing \n'Respect for Elders'")+
  ylim(0,1)+ theme_bw()

plot_respect

Respect.aov<-aov(Respect~as.factor(treatment), data=limited)
summary(Respect.aov)
TukeyHSD(Respect.aov)
eta_squared(Respect.aov, partial = FALSE)



obed_mean <- summarySE(limited, measurevar="Obedient", groupvars="treatment2",
                       na.rm=TRUE, conf.interval=.95)

obed_mean
obed_mean$test<-c("Control", "Boy", "Girl")

plot_obed<-ggplot(obed_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl', 'After', 'Before')), y=Obedient, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Obedient-ci, ymax=Obedient+ci), colour="red",data=obed_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Obedient'")+
  ylim(0,1)+ theme_bw()

plot_obed

Obedience.aov<-aov(Obedient~as.factor(treatment), data=limited)
summary(Obedience.aov)
TukeyHSD(Obedience.aov)
eta_squared(Obedience.aov, partial = FALSE)




GM_mean <- summarySE(limited, measurevar="GoodMannered", groupvars="treatment2",
                     na.rm=TRUE, conf.interval=.95)

GM_mean
GM_mean$test<-c("Control", "Boy", "Girl")

plot_GM<-
  ggplot(GM_mean, aes(x=factor(test, 
                               level =c('Control','Boy', 'Girl', 'After', 'Before')), 
                      y=                        GoodMannered, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, 
                aes(ymin=GoodMannered-ci, 
                    ymax=GoodMannered+ci), 
                colour="red",data=GM_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Good Mannered'")+
  ylim(0,1)+ theme_bw()

plot_GM

GM.aov<-aov(GoodMannered~as.factor(treatment), data=limited)
summary(GM.aov)
TukeyHSD(GM.aov)




WB_mean <- summarySE(limited, measurevar="WellBehaved", groupvars="treatment2",
                     na.rm=TRUE, conf.interval=.95)

WB_mean
WB_mean$test<-c("Control", "Boy", "Girl")


plot_WB<-ggplot(WB_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), 
                             y=WellBehaved, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=WellBehaved-ci, 
                              ymax=WellBehaved+ci), colour="red", 
                data=WB_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Well Behaved'")+ylim(0,1)+ theme_bw()
plot_WB+ theme_bw()

WB.aov<-aov(WellBehaved~as.factor(treatment), data=limited)
summary(WB.aov)
TukeyHSD(WB.aov)

Polite_mean <- summarySE(limited, measurevar="Polite", groupvars="treatment2",
                         na.rm=TRUE, conf.interval=.95)

Polite_mean
Polite_mean$test<-c("Control", "Boy", "Girl")


plot_Polite<-ggplot(Polite_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), 
                                     y=Polite, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Polite-ci, 
                              ymax=Polite+ci), colour="red", 
                data=Polite_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Polite'")+ylim(0,1)+ theme_bw()
plot_Polite+ theme_bw()

Polite.aov<-aov(Polite~as.factor(treatment), data=limited)
summary(Polite.aov)
TukeyHSD(Polite.aov)

ORD_mean <- summarySE(limited, measurevar="Orderly", groupvars="treatment2",
                      na.rm=TRUE, conf.interval=.95)

ORD_mean
ORD_mean$test<-c("Control", "Boy", "Girl")


plot_ORD<-ggplot(ORD_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), 
                               y=Orderly, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Orderly-ci, 
                              ymax=Orderly+ci), colour="red", 
                data=ORD_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Orderly'")+ylim(0,1)+ theme_bw()
plot_ORD

ORD.aov<-aov(Orderly~as.factor(treatment), data=limited)
summary(ORD.aov)
TukeyHSD(ORD.aov)

DIS_mean <- summarySE(limited, measurevar="Disciplined", groupvars="treatment2",
                      na.rm=TRUE, conf.interval=.95)

DIS_mean
DIS_mean$test<-c("Control", "Boy", "Girl")


plot_DIS<-ggplot(DIS_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), 
                               y=Disciplined, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Disciplined-ci, 
                              ymax=Disciplined+ci), colour="red", 
                data=DIS_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Disciplined'")+ylim(0,1)+ theme_bw()
plot_DIS

DIS.aov<-aov(Disciplined~as.factor(treatment), data=limited)
summary(DIS.aov)
TukeyHSD(DIS.aov)

Loyal_mean <- summarySE(limited, measurevar="Loyal", groupvars="treatment2",
                        na.rm=TRUE, conf.interval=.95)

Loyal_mean
Loyal_mean$test<-c("Control", "Boy", "Girl")


plot_Loyal<-ggplot(Loyal_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), 
                                   y=Loyal, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Loyal-ci, 
                              ymax=Loyal+ci), colour="red", 
                data=Loyal_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Proportion choosing \n'Loyal'")+ylim(0,1)+ theme_bw()
plot_Loyal

Loyal.aov<-aov(Loyal~as.factor(treatment), data=limited)
summary(Loyal.aov)
TukeyHSD(Loyal.aov)



library(cowplot)
pdf(file="Figure3a.pdf")

plot_grid(plot_respect, plot_obed, plot_WB, plot_GM)
dev.off()
pdf(file="Figure3b.pdf")

plot_grid(plot_Polite, plot_ORD, plot_DIS, plot_Loyal)
dev.off()

##MAKE THE AUTHORITARIANISM MEASURE

study3$Author<-study3$Respect+study3$Obedient+
  study3$GoodMannered+study3$WellBehaved+study3$Polite+
  study3$Orderly+study3$Disciplined+study3$Loyal

limited$Author<-limited$Respect+limited$Obedient+
  limited$GoodMannered+limited$WellBehaved+limited$Polite+
  limited$Orderly+limited$Disciplined+limited$Loyal

Author_mean <- summarySE(limited, measurevar="Author", groupvars="treatment2",
                         na.rm=TRUE, conf.interval=.95)

Author_mean
Author_mean$test<-c("Control", "Boy", "Girl")


plot_Author<-ggplot(Author_mean, aes(x=factor(test, level =c('Control','Boy', 'Girl')), 
                                     y=Author, group=1)) +
  geom_line() +
  geom_errorbar(width=.1, aes(ymin=Author-ci, 
                              ymax=Author+ci), colour="red", 
                data=Author_mean) +
  geom_point(shape=21, size=1, fill="white") +
  xlab("Treatment")+
  ylab("Authoritarianism measure")+ylim(0,5)+ theme_bw()

pdf(file="Figure6.pdf")

plot_Author
dev.off()
Author.aov<-aov(Author~as.factor(treatment), data=limited)
summary(Author.aov)
TukeyHSD(Author.aov)



library(cowplot)
plot_grid(plot_respect, plot_obed, plot_WB, plot_GM)

plot_grid(plot_Polite, plot_ORD, plot_DIS, plot_Loyal)



#plot_grid(plot_auth, plot_auth_c, plot_auth_b, plot_auth_g,
#          labels = c('Full Sample', 'Control', 'Boy', 'Girl'))



